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^ . Abstract 
O 

We study a bound system of N identical bosons interacting by model pair po- 
Q-r tentials of the form V(r) = A sgn(p)r p + fi, A > 0, ^ > 0. By using a variational trial 

' function and the 'equivalent 2-body method', we find explicit upper and lower bound 

formulas for the iV-particle ground-state energy in arbitrary spatial dimensions d > 3 
for the two cases p = 2 and p = -1. It is demonstrated that the upper bound can 
be systematically improved with the aid of a special large-TV limit in collective field 
^ ! theory. 

PACS 03.65.Ge; 31.15.Pf; 03.75.Fi. 
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1. Introduction 

The principal subject of this paper is the ground-state energy of a system com- 
posed of N identical bosons bound together by pair potentials. Such systems usually 
collapse in the large- N limit; that is to say, the binding energy per particle rises with 
N to infinity [1-2.] In two earlier papers we studied gravitating boson systems [3], and 
systems of bosons or fermions which interact by wider classes of purely attractive pair 
potential [4]. In the present paper we extend this work for the boson case to systems 
with pair potentials which bind the system but have a soft repulsive core of the form 
^. The methods we develop to analyse these systems may have application to models 
for Bose-Einstein condensates [5, 6] in which there is at present much renewed inter- 
est [7-17]. As in our earlier work, the goal of the present study is to provide simple 
general energy bounds as functions of the parameters m, A, fj,, and N. Since is not 
difficult to do so, we also allow for an arbitrary number d > 3 of spatial dimensions. 

Formulas for TV-particle upper and lower energy bounds are derived for pair po- 
tentials Vij = V a f(rij) whose shapes /(r) are either harmonic oscillators with a soft 
core f(r) = Xr 2 + -j^, or Kratzer potentials f(r) = — £ + ^. The shapes of these two 
potentials are shown in Fig.(l) for A = \i = 1. The first of the models could perhaps be 
considered as a generalization of the soluble 1-dimensional Calogero model [18- 20] to 
dimension d > 3. For such potentials, with minima /(f) at r = f > 0, provided there 
are enough spatial dimensions available (d+ 1 > N) for every pair distance to satisfy 
nj = f, then, in classical mechanics the lowest energy would be equal to (^)V f(r). 
This expression provides a lower energy bound in both classical and quantum me- 
chanics; our general energy bounds (valid for all N > 2 ) will show that this value is 
approached asymptotically for large V a in the limit d — » oo; this is possible because 
the positive contribution from zero-point oscillations varies like (N - l)(NV )% and is 
eventually dominated by the static potential term. 

The non-individuality of identical quantum-mechanical particles introduces a very 
powerful constraint which allows us to relate the iV-body problem to a specially con- 
structed (reduced) 2-body problem having an overall factor of N - 1 and a potential 
enhanced by N/2, corresponding respectively to the N — l relative kinetic-energy op- 
erators and the N(N — l)/2 pair potentials. This 'equivalent two-body' notion, which 
is central to our approach to the study of these systems, will be formulated explicitly 
in Section (2). The TV-body energy £ is described by function E = F N (v) = £/(N - 1), 
where v = %V , and V a is the coupling parameter (in units K = m = l). Since v 
includes the factor N, it follows that if we consider & finite value of v and let TV — > oo, 
this implies that V a — » 0. The well-known exactly-soluble ^-potential problem in 
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one dimension [21] provides a convenient illustration of the energy function F(v). If 
V(r) = —V 5(r), we have exactly 



1 

N 



v\ N>2, v=^V a . (1.1) 



The shape of the potential dictates the form v 2 of these energy curves, and there is 
a distinct curve for each value of N > 2. The curves of this particular example satisfy 
the general functional inequalities [3] 

F 2 {v) <E = F N {v) < F^v) < F^v) < F G (v), (1.2) 

where, as will be explained later, F^v) is an upper bound obtained with the aid 
of collective- field theory in the large- N limit, and F G (v) is the upper-bound curve 
obtained by employing an TV-body translation-invariant Gaussian trial wave function. 
The ordering of the F N (v) curves with TV is a consequence of the monotonic increase 
in the severity of the boson-symmetrization constraint with increasing N. For the 
soluble J-potential model we have [3] explicitly: F^l) = -1/6, F (l) = -0.164868, and 
F G (l) = -1/27T. The family of densities used for the variational collective-field upper 
bound is given by 

«/>(r) = e~^ q , b, q>0. (1.3) 

The energy is minimized with respect to the positive scale and power parameters b 
and q. The effectiveness of this useful 2-parameter family of trial densities is also 
demonstrated by the applications discussed in Sections (4) and (5) below. We con- 
jecture that F^v) is close to F^v) for the more general problems considered in this 
paper but we know no way of proving such a claim at this time. For another well- 
known soluble problem, the d- dimensional harmonic oscillator V(r) = V a r 2 , with , 
d > 1, the inequalities (1.2) all collapse together to the common exact value E = dvi. 

2. The iV-body problem in the centre-of-mass frame 

The Hamiltonian, with center-of-mass removed, for a system of TV identical par- 
ticles each of mass m interacting via central pair potentials may be written 

AT / N \ 2 N / 1 ,s 

i=l V=l / 3>i=l V 7 

where V Q and a are respectively the depth and range parameters of the potential with 
shape /. By algebraic rearrangement (2.1) may be rewritten in the more symmetrical 
form 

«- t {^<p.-p-) 2 +^(^)}- 

j>i=l k \ / J 
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We now define new coordinates by p = BR, where p = [pi] and R = [r,] are column 
vectors of the new and old coordinates, respectively, and B is a real constant NxN 
matrix. For convenience we require all the rows of B to be unit vectors, we let the 
elements of the first row all be equal to so that pi is proportional to the centre- 
of-mass coordinate; we also require that the remaining N—l rows of B be orthogonal 
to the first row, so that they define a set of N - 1 relative coordinates. One more row 
is also fixed so that we have at least one pair distance at our disposal, namely 

For boson systems, we have found that Jacobi relative coordinates, for which B is 
orthogonal, are the most useful. Thus, corresponding to the transformation p = BR 
of the coordinates, it follows that the column vector P of the associated momenta 
transforms to the new momenta n = [7^] by the relation n = (B T y 1 P = P. If * is 
any translation-invariant wave function for the TV-body system composed of identical 
bosons, then we can write [22, 23] the following mean energy relation between the 
iV-body and 2-body systems: 

(*,W*) = (*,H*), (2.4) 
where the 'reduced' two-particle Hamiltonian IH is given by 

Further simplifications can be achieved if we work with dimensionless quantities. 
We suppose that the translation-invariant TV-body energy is £ and we define the 
dimensionless energy and coupling parameters E and v by the expressions 

m£a 2 NmV a 2 

E=- — v= 5 — . (2.6) 

(N-l)h 2 2h 2 

It is then natural to define a dimensionless versions of the reduced 2-body Hamiltonian 
IH and the relative coordinate p 2 by the relations 

H= ( ^ H ^ 2 =-A + »/(r), r = V2p 2 /a=(r 1 -r 2 )/a, r=||r||. (2.7) 

We note that the Hamiltonian H depends on N only through the dimensionless 
coupling parameter v. By the Rayleigh-Ritz (min-max) principle [24-26], we have the 
following characterization of the TV-body ground-state energy parameter E in terms 
of H -. 

Pity) 

E = ^-^=F N (v), (2.8) 
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where * is a translation-invariant function of the JV - 1 relative coordinates (and 
spin variables, if any) which is symmetric under the permutation of the JV individual- 
particle indices. The TV-body energy £ is recovered from E by inverting (2.6). Thus 
we have explicitly: 



However, for the remainder of this paper we shall work with the dimensionless form 
E = F N {v) -. the problem is to find or approximate F N . In each situation we shall have 
to specify the dimension d of physical space; usually d = 3. . 

3. Energy bounds 

The energy bounds used in this paper are summarized in terms of the F functions 
by Eq.(1.2). The history of the equivalent 2-body method for boson systems has been 
described in the earlier papers [3, 4] and in the references therein. The main result is a 
general energy lower bound which, for boson systems with orthogonal Jacobi relative 
coordinates, is given by 



where F 2 (v) is the lowest eigenvalue of the 1-particle ('reduced' 2-particle) Hamilto- 
nian H = -A+vf(r). With equal simplicity, our weakest upper bound F G (N), provided 
by a Gaussian boson trial function may also be expressed in terms of the 1-body 
operator H. This is a consequence of the following argument. If and only if the sym- 
metric translation-invariant function * is Gaussian [27, 28], it may be factored in 
the form 



But the equivalence (2.4) then implies, in this case, that E < (ip,Hip)\\tp\\^ 2 . This 
explains why the inequalities (1.2) collapse together in the case of the harmonic 
oscillator for which the exact l-body lowest eigenfunction of H is also Gaussian. In 
this argument we assume for the upper bound that <H> has been optimized with 
respect to the scale of the wave function. The boson symmetry of these Gaussian 
functions is demonstrated most clearly by the following algebraic identity: 




(2.9) 



F 2 (v) <E = F N (v) 



(3.1) 



*(P2,P3, ...,Pn) = ^(p 2 )v{P3,P4, ■ ■ -,Pn)- 



(3.2) 



N 



N 



S ( r * - r j 



) 2 = n^pI 



(3.3) 



j>i=l 



The potentially better upper bound F 4> (N) is more complicated both to derive 
and to compute [3] . It is found by considering the collective field model in the large- JV 
limit. The energy so obtained is an upper bound to F^v) and this, in turn, may be 
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closely approximated from above by optimizing the right-hand side of the following 
equation: 

F N (v) < F^v) = l -J {V f^ 2 d d r + vj J 0(r)/(|r - r'\)cj>(r')d d rd d r l , (3.4) 

where <f>(r) is a trial probability density function (for inter-particle distances) satis- 
fying / 4>(r)d d r = 1. It has been shown [3] that if <j> is Gaussian, and the energy is 
optimized with respect to a scale parameter, then the result F^{v) is identical to the 
upper bound F G (v) obtained with the aid of a scale-optimized translation-invariant 
Gaussian trial wave function. The wider possible choice of the form of the probability 
density 4> allows us to transcend the Gaussian bound whilst still working essentially 
with a variational 'function' of single variable for the iV-body problem. 

In order to facilitate the reproduction of our results and the application of the 
method to other problems, we make the following explicit technical remarks concern- 
ing the case d = 3. It is very helpful to think of the probability density as a function 
(f>(r) = w(r/b) = w(s), which depends on the remaining parameter q, to be discussed 
later. If we let I = / °° w(s)s 2 ds, then the kinetic energy integral becomes 

<KE > * fMMl!^ (3.5) 

8/fo 2 7 n w(s) v ; 

and potential energy integrals for pure powers may be written 

up r°° r°° 

<r p >= W( \ dtw(t)t d SW ( S ) S {(s + tr+ 2 -(s-tT+ 2 }, pjt-2. (3.6) 



t 



I 2 (p + 2) J 
For the case p = -2 we have instead 



1 1 f°° f°° ( s I t\ 

< ^2 > = j2tfj dtw(t)tj dsw(s)sln( — -J . (3.7) 

With the terms expressed in this form, the minimization over scale b > can often 
be carried out explicitly yielding an algebraic expression which then needs to be 
minimized with respect to the remaining parameter q. We have explored various 
alternative forms for w(s), such as s 9 exp(-s 2 ) and exp(-(s - q) 2 ), but have found it 
difficult to improve on the variable-power family w(s) = cxp(-s q ) which we have used 
to obtain the results discussed in detail in Sections (4) and (5) below. 

4. The harmonic oscillator with a soft core 

We now choose the potential shape to be 

/(r) = Ar 2 + ^, A>0, M > 0. (4.1) 
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The lower bound F 2 (v) is given by the lowest eigenvalue of H — -A + vf(r) in d 
dimensions, a problem which is discussed, for example, in Ref. [29]. The Gaussian up- 
per bound F G (v) is provided by minimizing the Rayleigh quotient (ip,Hi/>)/(il),il>) with 
respect to the scale variable a in ip(r) = e~i ar ' 2 in d dimensions. These calculations 
yield the following bounds on the TV-boson energy parameter E (the energy itself is 
recovered essentially by multiplying E by N — 1, according to Eq.(2.10)) 



< E 



< 



(dv\y- 



4v/i 
d-2 



d>3. 



(4.2) 



The Gaussian trial wave function allows us to compute an approximate value for the 
mean-squared pair separation a 2 , a measure of the size of the system. We find 



<(ri - r 2 ) 2 > 



d 



2(v\y 



i + 



4v/j, 



d{d-2) 



(4.3) 



We now look at some special cases. For the harmonic oscillator, /x = 0, the in- 
equalities collapse to the exact value E = d(vX)i. The asymptotic forms of the bounds 
as => v — > oo and \i > are given by 



2v(A^) 3 < E < - 2v(Xn)i 



d-2 



(4.4) 



Thus for large TV the energy per particle E increases like N, although the 'size' of the 
system (as estimated by the Gaussian wave function) approaches the constant value 
a = (d/(d- 2))i(/i/A)3. If d is now taken large, the asymptotic form of the energy 
approaches the classical expression E = vf(f), where f = a = (/x/A)i is the position of 
the minimum of f(r). 

For the special case d = 3 we obtain in general 



(uA)2 2 + {l + 4fiv) 2 < E < («A)2 [9 + 12v/j] z 



(4.5) 



These results are shown in Fig. (2) along with the improved upper bound ('dashed' 
curve) obtained with the trial probability density function (1.3). The extreme optimal 
q(v) values in the range shown were: q{2) = 2.8593, g(20) = 4.460. Hence the optimal 
probability density <f>{r) is found to be quite far from Gaussian. 

5. The Kratzer potential 

The Kratzer potential has shape function given [30-39] by 

/(r) = + A>0, M >0. (5.1) 
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The iV-particle lower bound F 2 (v) is provided by the lowest eigenvalue of H = 
—A + vf(r) in d dimensions, a problem which is discussed, for example, in Ref. [40]. 
Meanwhile an exactly similar calculation to that described in Section (4), with a 
Gaussian trial function, generates the corresponding upper bound. We find in this 
way that the energy parameter E satisfies the following inequalities: 



(vX) 2 



1 + {(d- 2) 2 +4vnY- 



< E < — 



2d 



8v^ 
d-2 



d>3, 



(5.2) 



where the d-dependent constant -y d is given by 



id 



r(f) ' 



73 



2 

7*' 



(5.3) 



The Gaussian trial wave function again allows us to compute an approximate value 
for the mean-squared pair separation. We find 



<(ri - r 2 ) 2 > 



2(vX ld ) 2 \~ ' d(d-2) / 
For the pure gravitational case n = the energy inequalities become 



(vXf 



< E < 



(d-iy 

The asymptotic forms of the bounds as 

~-t < E < - 



2d 

oo and n > are given by 
where M(d) - \ ^( d ~ 2 ) 



(5.4) 



(5.5) 



(5.6) 



The function M(d) increases monotonically with d to 1; M(3) = 2/tt w 0.63662; M(8) > 
0.9. Thus, as for the previous model, the energy per particle E increases like N. 
Meanwhile the size (as estimated by the Gaussian wave function) approaches the 
constant value a = (2V2d/(^ d (d - 2))(ji/\). If d is now taken large, the asymptotic 
form of the energy (again) approaches the classical expression E = vf(f), in which 
f = a = 2/j,/\ is the position of the minimum of /(r). 
For the special case d = 3 we obtain in general 

(vX) 2 „ ^ ^ (vXf 



1 + (1 + 4v/j,) 2 



< E < 



7T [| + 2VH] ' 



(5.7) 



The bounds for three dimensions are shown in Fig. (3), along with the improved upper 
bound ('dashed' curve) obtained with the trial probability density function (1.3). 
The extreme optimal q(v) values in the range shown were for this problem: q{2) = 
2.0017, q(20) = 3.237. The optimal probability density (in this family) is found to be 
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almost Gaussian for small v < 3 but very different from Gaussian for larger values of 
v. 

6. Conclusion 

The main purpose of this paper is to derive general bounds for the energy of N- 
boson systems which are bound together by pair potentials with soft cores. We have 
examined two such models and we have provided upper and lower bound formulas 
valid for all N > 2 and d > 3, and for all values of the potential parameters that 
bind the system. These bounds are expressed in terms of the energy parameter E = 
£/(N - 1) = F N (v), where (in units with h = m = l) v = NV Q /2. If the potential shape 
is such that F 2 (v) is close to F^v), then we obtain close upper and lower bounds 
valid for all N > 2. The upper bound F G (v) provided by a Gaussian trial function 
may be improved to F^v) which is derived by using a (possibly non-Gaussian) trial 
probability density 4> in a limiting form of collective field theory. In order to do 
significantly better one might hand craft a trial wave function for a particular f(r) 
and N. Of course, to be secure about this wave function one would still need to find 
a good lower bound, a goal usually not easy to achieve. The global results we have 
obtained show that these boson systems are asymptotically bounded by expressions 
of the form £ — cN 2 ; if d is also large, we have shown that the heuristic classical 
expression £ = (^)V f(r) is reached asymptotically. Since the many-body problem 
continues to offer a serious challenge for direct numerical solution, it is very helpful 
to have some explicit analytic upper and lower energy bounds and size estimates for 
model systems such those discussed in this paper. 
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Figure 1. 

The shapes of the two potentials studied: the harmonic oscillator with soft core, 
H = f(r) = Ar 2 +^r~ 2 , and Kratzer's potential, K = f(r) = -Ar _1 +,ur _2 , with A = fj, = 1. 
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Figure 2. 

Bounds for the energy parameter E(v) = £/(N - 1) of the iV-boson system with pair 
potentials of the form /(r) = V„(r 2 +r~ 2 ), as a function of v = NV /2. The graphs show: 
U the upper bound F G (v) found by a Gaussian trial function, L the lower bound F 2 (v) 
by the equivalent 2-body method, and 'dashed' the upper bound F^{v) by collective 
field theory. 



General energy bounds for systems of bosons with soft cores 



page 14 



E(v) 




2 4 6 8 10 12 14 16 18 20 



Figure 3. 

Bounds for the energy parameter E(v) = £/{N - 1) of the iV-boson system with pair 
potentials of the form /(r) = V^-r -1 + r~ 2 ), as a function of v = NVq/2. The graphs 
show: U the upper bound F G (v) found by a Gaussian trial function, L the lower 
bound F 2 (v) by the equivalent 2-body method, and 'dashed' the upper bound F^v) 
by collective field theory. 



